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ABSTRACT 

Dual-Frequency Polarized Scatterometer (DFPSCAT) is a pencil-beam rotating 
scatterometer which is designed for snow water equivalent (SWE) measurement, and 
Doppler Beam Sharpening (DBS) technique is proposed for DFPSCAT to achieve the 
azimuth resolution. However, the DBS technique is inapplicable for the 
forward-looking and afterward-looking region. Based on an approximate aperiodic 
model of scatterometer echo signal, an improved adaptive regularization 
deconvolution algorithm with gradient histogram preservation (GHP) constraint is 
implemented to settle the problem. To investigate its performance of resolution 
enhancement and resulted accuracy, both a synthetic backscattering coefficient (a°) 
field reconstruction and SWE o° reconstruction are carried out. The results show that 
the proposed method can recover the truth signal and achieve azimuth resolution of 
2km with the designed scatterometer system, which is required by the SWE retrieval. 
Moreover, the relative errors of reconstructed o° are less than 0.5dB that satisfy the 
accuracy requirement for SWE retrieval, and comparisons with observed results 
shows that the error reduction are more than 0.03dB. Meanwhile, a comparison 
between the proposed algorithm and some existing resolution enhancement methods 
is analyzed, which concludes that the proposed method can obtain a comparable 
resolution enhancement as Lı method and has less noise. The technique is also 
verified with ASCAT scatterometer data. 
Index Terms: Dual-Frequency Polarized Scatterometer (DFPSCAT), snow water 


equivalent (SWE), azimuth super-resolution, reconstruction accuracy 


I. INTRODUCTION 

S paceborne radar scatterometers mainly provide wide-area measurements of 

ocean surface wind vectors and have been widely used in the investigations of 
oceanography and meteorology. With advantages of wide swath, quick global 
coverage and high measurement accuracy, scatterometers have become the desirable 
candidates for land and ice applications [1]. For the conventional scatterometers, 
resolution in range dimension is processed into a few kilometers using the range gate 
and that in azimuth dimension is still determined by the antenna beamwidth about a 
few tens kilometers. For land and ice applications, spatial resolution comparable with 
visible/IR imaging radiometers, i.e., the 1-5km resolution Advanced Very High 
Resolution Radiometer (AVHRR), is demanded [2]. 

To achieve the required resolution, meanwhile keep the advantage of wide swath, 
several high resolution reconstruction methods have been proposed for scatterometer 
image reconstruction, such as additive algebraic reconstruction technique (AART), 
multiplicative algebraic reconstruction technique (MART) and scatterometer image 
reconstruction (SIR) [3, 4]. As an improved version of MART, SIR has the fastest 
convergence and the least noise amplification, and has been applied to the 
scatterometric data for some specific applications such as ice mapping and iceberg 
tracking et al. In addition, Williams established a Bayesian-based scatterometer 
reconstruction framework, in which the scatterometer noise was taken into 


consideration [5]. The maximum posteriori (MAP) estimator which is developed by 


Williams can obtain consistent results as the well-established SIR algorithm, but 
significantly enhances the resolution at the expense of noise amplification. 

The Water Cycle Observation Mission (WCOM) is proposed to improve the 
capability of synergetic observation of key water cycle parameters. Dual-Frequency 
Polarized Scatterometer (DFPSCAT) is one of the payloads, proposed for WCOM 
along with two others, the L-S-C tri-frequency Interferometeric Microwave Imager 
(IMI) and the Polarimetric Microwave Imager (PMI). DFPSCAT is dedicated to 
retrieve snow water equivalent (SWE) by simultaneously measuring surface and 
volume scattering from snow coverage [6, 7]. According to the requirements of SWE 
retrieval, the spatial resolution of DFPSCAT should be less than 2km to mitigate the 
retrieval errors related to the mix-pixel problem, and its radiometric precision should 
be better than 0.5dB [8]. Since the existing scatterometer systems and technology 
cannot satisfy the resolution and accuracy demands, Doppler Beam Sharpening (DBS) 
technique is proposed for DFPSCAT to achieve the goal of the resolution in azimuth 
direction. However, the DBS technique is inapplicable for the forward-looking and 
afterward-looking region, where the Doppler discrimination is confused with range 
discrimination [9]. Therefore, it is necessary to find some other super-resolution 
methods to overcome the difficulty. 

Guan et al. derived a forward-looking radar echo model for airborne radar and 
concluded the radar echo was equivalent to a convolution operation between target 
normalized radar backscatter coefficient (0°) and antenna pattern [10]. Then, the echo 


model can be solved by some deconvolution methods to achieve high resolution. With 


this motivation, the sampling rate of DFPSCAT in azimuth direction is increased to 
meet the desired resolution requirement, and then the deconvolution methods are 
utilized to improve the azimuth resolution for the forward- and after-looking region. 
However, the noise and the nulls of antenna pattern make the deconvolution become 
an ill-posed problem. Wiener filter eliminates the antenna nulls in frequency domain 
by introducing a constant regularization term, but the algorithm requires high 
signal-to-noise ratio (SNR) [11]. To solve the problem, Richards proposed an 
constrained iterative algorithm to reduce the demand of SNR, but it is still strongly 
influenced by noise [12]. Guan introduced a statistical optimization method to reduce 
the noise sensitivity, but the method amplified noise locally during the iteration [10]. 
Zou proposed a norm regularization method used in radar azimuth super-resolution, 
and obtained promising results [13]. Both the L and Lı norm regularization model 
were set up, and the results turn out that the latter one could obtain more effective 
resolution improvement especially for sparse target signal. Besides, Wang utilized a 
total variation regularization (TV) method for scatterometer image reconstruction of 
the rotating fan-beam scatterometer (RFSCAT), and achieved an improved quality of 
resolution enhancement and noise reduction [14]. 

The aim of this paper is to develop a theoretical framework that allows the 
scatterometer forward-looking echo signal model to be appropriately solved by 
deconvolution method. Two kinds of approximate model (circulant matrix model and 
aperiodic matrix model) are introduced. Considering the over smoothness of the 


adaptive regularization method, an improved adaptive regularization deconvolution 


algorithm with gradient histogram preservation (GHP) constraint has been 
implemented to solve the echo model. The improved adaptive regularization 
algorithm integrates gradient histogram preservation into the adaptive regularization 
method which can better perverse the fine variations of o°. To illustrate the feasibility 
of the deconvolution framework and the effectiveness of the proposed algorithm, both 
simulation and actual scatterometer data validation are conducted. The performance 
of resolution enhancement and processing accuracy of the proposed method is 
investigated, and compared with that of conventional approaches. 

This paper is organized as follows. Section II deduces a DFPSCAT echo signal 
model for super resolution processing. Section III proposes an improved adaptive 
regularization deconvolution algorithm and describes its iterative specification. 
Section IV discusses the performance of resolution enhancement and processing 
accuracy for the proposed algorithm. Section V validates the feasibility of the 
deconvolution framework and the effectiveness of the proposed algorithm using real 
scatterometer data. Finally, Section VI concludes this paper. 

II. DFPSCAT ECHO SIGNAL MODEL 
A. THE INSTRUMENT 
DFPSCAT is a real aperture radar operating at X-band (9.6 GHz) frequency and 
Ku-band (17 GHz) frequency using both vertical and horizontal polarization (i.e., VV, 
HH and HV polarization). It transmits long pulses with linear frequency modulation, 
and the received echoes are processed with a de-chirping technique to achieve the 


range resolution enhancement. The original range resolution is about several-hundred 
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meters, with transmitted signal bandwidth of about 5MHz, and can be further 
averaged into 2-5km spaced slices along-elevation before being downlinked to ground. 
The DBS technique is utilized in the outer swath regions of DFPSCAT, and results in 
the azimuth resolution of 1-3km. For the forward- and afterward-looking region of 
200-300km, the azimuth resolution is still determined by the antenna beamwidth in 
the azimuth direction. Fig. 1 illustrates the observation geometry of DFPSCAT. The 
center of each pencil beam has a nominal off nadir angle of 39° together with the 
altitude of 600 km to ensure that the swath is about 1000km. In the figure, “High 
Resolution” denotes the DBS processing region, and “Low Resolution” represents the 


forward- and afterward-looking region for the super resolution algorithm processing. 


Nadir Track 


Antenna Beam 


sa gg 


C] Low Resolution 


LJ High Resolution 
Nadir gap in high-resolution 
radar data:200-300km 


Fig. 1 The observation geometry of DFPSCAT 


The preliminary designed system parameters of DFPSCAT are shown in Table 1. 
The antenna 3dB beamwidths are 1.2° and 1° for X-band and Ku-band, respectively. 
Considering the accuracy of SWE retrieval which is greatly influenced by spatial 


resolution, the resolution of 2-5km for DFPSCAT is required. Based on the previous 


illustration, the period repeat frequency (PRF) of DFPSCAT is designed to satisfy 


Ikm sampling rate in azimuth direction. 


Table 1 DFPSCAT system parameters 


Parameters X-band Ku-band 
Frequency 9.6 GHz 17GHz 
Polarization VV, HV 
Altitude 600km 
Incidence Angle 43.5° 
PRF 1320Hz 
Rotation Rate 15.4 revolutions per minute (rpm) 

Swath 1000km 

3dB Beamwidths 1.2°x1.2° 1.0°x1.0° 


3dB Footprint Dimensions 
22.5x16.4km 16.4x14km 
(elevationxazimuth) 


B. ECHO MODEL 
From Fig. 1, we can see that the echo of DFPSCAT varies with slant range R and 


azimuth angle ø. The amplitude of the echo is the result of filtering and weighting of 


the backscattering distribution of the observed targets. Assuming the transmit signal is 


a rectangular pulse signal 


p(t) =rect(r) ef C7 (1) 
where t represents transmit time, rect(z) denotes the rectangular pulse envelope, feis 
carrier frequency and ¢(z) is phase factor. 


The echo signal after pulse compression can be expressed as 


Z gR 
y(R, p) = Co (R, p) ® HR pač] @) sine 28 a” a # c i 
DP C C 


where C, is the radar system constant, o° is the normalized radar backscatter 
coefficient, which is the function of the slant range R and the azimuth angle g, h 
represents antenna pattern, A is the electromagnetic wave length, c is the light speed, 
® represents convolution symbol. 

On the forward-looking region of DFPSCAT, the slant ranges are the same for 
the same slice. Then, Eq. (2) can be simplified as 

WP) =(P) Dh) (3) 
Fig. 2 shows the azimuth sampling schematic of DFPSCAT. The continuous 


echo signal y can be discretized as 
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where N and M indicate the sample numbers of y and A, respectively. 
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Fig. 2 Azimuth sampling schematic of DFPSCAT (R = Ro) 


Due to the number of unknown o° exceeds the number of equations, there is no 
definite solution of Eq. (4). To settle the problem, one way is to utilize a circulant 
approximate model, which first spares the front M-1 columns of the matrix in Eq. (4), 
and then modifies the remaining part into a periodic structure[15]. The corresponding 


circulant model is 


[ h hy -1 7 h, | 
h h h, Gy 
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The other way is to modify the partial convolution of Eq. (4) into a complete 
convolution[15], namely the aperiodic approximate model. Then, Eq. (4) can be 


modified as 


hy 
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Compared with the circulant model, the aperiodic model is overdetermined 

problem, which is more favorable for solving the problem. Therefore, we adopt the 
aperiodic model and the matrix-vector form which can be expressed by 

y, =Ho” (7) 

where H and o” denote the matrix and vector in Eq. (6), respectively. Thus, y, can be 


obtained through the following expression 


Ja =Y, =H i| o ya a o” O as 0 Cr 1 Ni (8) 
where Ho represents the convolution matrix of Eq. (4). In fact, the boundary value of 
olua and Oyma Oya are unknown in real situation. Tekalp proposed to 
use the observed values to replace the boundary values, which produced large errors 
at the boundary of y, [16]. Instead, the paper proposes to utilize the previously 
described circulant model to calculate the unknown values in Eq. (8). The results 
show that y, obtained through our proposed method is closer to the complete 
convolution results. 

For scatterometry, the measurement error mainly consists of the following three 
components: transmission error caused by the fading effect of radar and thermal noise 
of receiver, calibration error resulting from the uncertainty of radar system parameters, 
and model error introduced since the inverse model and calibration model are 


inconsistent to the real scene. Usually, the value of the error is denoted by its 


normalized standard deviation, which is defined in terms of a so-called K, parameter. 


K = varLy, ] (9) 


aly] 


p 
For the scatterometer, the noisy measurement can be represented as 
y=(1+K,:N(0,))H o’ (10) 
where N(0,1) denotes the standard normal distribution, K, denotes the normalized 
standard deviation of measurement error. 
M. REGULARIZATION ALGORITHM 
Based on the former analysis, the azimuth super resolution is equivalent to the 


deconvolution problem. Due to the ill-poseness of deconvolution, one popular solving 


method is the regularization-based technique which makes the equation solutions 
continuously depend on observed data through introducing some regularization 


operators[17]. Then the desired o° could be obtained by 


6" =arg min| |Ho"—y n +a Y(o”)|,a>0 (11) 
where ¥ denotes some regularization term which depends on o° priors, œ is the 
regularization parameter. 
A. SELECTION OF Y 

In most practical applications, scatterometer image reconstruction is ill posed. 
Some form of regularization is required to estimate o° image from the scatterometer 


measurements. The proposed AART, MART, SIR and MAP algorithms all have a 


built-in regularization term, specifically 


AART: Plo. 

L 
MART(SIR): W=(0")' log, (o°) 
MAP: Y=] -0l 


These algorithms can be equivalent to the classical Tikhonov regularization 


algorithm (referred as TR algorithm) where ¥(o°)=|0° i 


in Eq. (11). Then, their 


e 
solution may be written as 
(H"H +a1)o° =H'y (12) 
Besides, an adaptive regularization method (denoted as AR) has been proposed 
to solve inverse problem of linear equations, and shows some promising results[18]. 
With the prerequisite that H’H is reversible, the specific form of the algorithm is 


((H"H) +al)o° =H"H H’y (13) 


Let s=H'y, ©=H'H, then the Eqs. (8), (12), (13) can be expressed in the 
form of Fø’ =s, then the corresponding F expression in Eqs. (8), (12), (13) can be 


2 
written as F=, F/*=O+alI and F“ ae i Comparing the above 


formulae, it is not difficult to find that when œ —> 0, both F and F” are 
approximate to F , whereas when ©® contains zero eigenvalue, F >0, F” >0 
and F” =a", which indicates Eq. (13) is a better approximation of Eq. (8). 
Therefore, the optimal solution of the echo model can be obtained by using the 
adaptive regularization algorithm. 

To evaluate the algorithm performance, one-dimensional (l-D) signal 
reconstruction simulation has been conducted, as shown in Fig. 3. The test signal is 
set as a narrow sinc function with a bandlimited spectrum. The antenna pattern is 
assumed as a wide rect function whose frequency response attenuates the high 


frequency components of the signal spectrum [3]. 


Amplitude 


Azimuth angle(deg) 


Fig. 3 Test signal and antenna pattern 
The reconstructed signal and its corresponding spectrum of different algorithms 


are investigated, as shown in Fig. 4. It is obvious that both AR and MAP algorithms 


can reconstruct the original signal with less bias compared with SIR and TR 


algorithms, whereas the AR algorithm has least noise amplification. 


SIR MAP SIR MAP 
1 1 0.015 0.015 
0.5 0.5 0.01 0.01 
0 0 0.005 0.005 
-0.5 -0.5 0 0 
-10 0 10 -10 0 10 0 500 0 500 
TR AR TR AR 
1 1 0.015 0.015 
0.5 0.5 0.01 0.01 
0 0 0.005 0.005 
0.5 -0.5 0 0 
-10 0 10 -10 0 10 0 500 0 500 


(a) Reconstructed signals (b) Spectrum of the reconstructed signals 


Fig. 4 Algorithm performance comparisons for noisy measurements 


B. AN IMPROVED ADAPTIVE REGULARIZATION METHOD 

Considering the above proposed L2-norm based adaptive regularization method might 
result in an over-smoothed solution, this section proposes an improved version of the 
adaptive regularization method. Since a good reconstructed o° should have a similar 
gradient distribution of the original ones, a gradient histogram preservation (GHP) 
term is integrated into the above adaptive regularization method [19]. The new 
proposed method is referred to as an improved adaptive regularization method with 
GHP constraint. 


The adaptive regularization method with GHP constraint can be expressed as 


A . 2 
6° =arg min |H o° -y| +a: 
oO 


o’ a 2 
Hopee] saans ‘8 
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where G denotes some odd function with monotonically non-descending property, V 


represents the gradient operator, hg is the histogram of la(vo") 


; hyo is the 


gradient histogram of o°, a and 2 are two regularization parameters. 


To solve Eq. (14), the gradient histogram of o°, h should be estimated 


Vo?” 


accurately. Zuo introduced an estimation method of reference gradient histogram 
from noisy image [19]. However, the method is only suitable for the denoising case. 
This paper substitutes o° with the gradient histogram of noise-free observation 


Ho’, and obtains the following adaptive regularization with GHP constraint model 


2 


2 
ata 


L 


eae |e» talets) -vali st. hg=h,, (15) 


where g=Ho’. 

C. GRADIENT HISTOGRAM ESTIMATION FOR Ho” 

To estimate the desired o° from y using Eq. (15), the reference histogram hy, 
should be determined, which is supposed to be the histogram of noise-free 
observation Vg in the paper. Substituting Ho’ with g, Eq. (10) can be written as 
y =g-+4n, where n denotes the measurement error and can be modeled as Gaussian 
white noise with variance 5°. Then, the relationship Vy =Vg+Vn_ can be obtained, 
where Vy are independent and identically distributed with the probability density 
function (PDF) of p,,. Since n subjects a Gaussian distribution of N (0, Oo), the 


corresponding PDF for Vn can be well approximate as follows[20]. 


1 (Vn) 
= -> 16 
Pyn 2 lS ex 452 ( ) 


Assuming Vg and Vn are independent, and Vg satisfies independent and 


identically distributed with PDF of Pyg, then py, can be expressed as 


Pvy (Vy =t) = | Py (V8 =4)* Pyn (Vn=t-a)da (17) 


Since the normalized histogram can be regarded as a discrete approximation of 
PDF, Py, Py,,and Py, canbe approximated by hy, hy, and hy, , respectively. 
Then Eq. (17) can be rewritten as 

Igy = hyg hyn (18) 
where © denotes the convolution operator. 


Using the above relationship, the estimation of hy 


g is equivalent to a 


deconvolution problem. Generally, the gradient histogram of Ho’ can be modeled 


as a generalized Gaussian distribution[21] with the following expression. 


_ y y 
Pv = aal (\Vg|/ 2) ) (19) 


where I represents complete Gamma function, / is a scale parameter, y is shape factor 
which is inversely proportional to the decreasing rate of the Gaussian distribution 
peak. Then the estimation of py, is converted into that of parameters [f, y]. 

D. ITERATIVE SPECIFICATION OF THE IMPROVED ADAPTIVE 

REGULARIZATION METHOD WITH GHP CONSTRAINT 

Considering the complexity of the optimization problem of Eq. (15), an iterative 
method based on steepest descent algorithm has been chosen. The iterative 
specification for the proposed mixed regularization method is summarized in 
Algorithm 1, which is described in Table 2. 


Define an object function as 


2 


0 
+9|G(Ve)—Vel|, (20) 
LC 


J =|Ho" yf, +a 


And its gradient can be written as 


T T 0 o” k ToT 
VJ =2H"(y-H'o"(k))+2a oA) «9H v” (G(Vg)-Vg) (21) 


T 


Then, the steepest descent algorithm satisfies the following relationship 
o° (k+1)=0°(k)+u-VJ (22) 
where y is the step size. 
Another problem is to select an optimal G(Vg) to satisfy the requirement of 
hg = hy . In the paper, G(Vg) is defined as follows [22]. 
G(Vg) =sgn(Vg)-T (|Vgl) (23) 
where T represents an feasible monotonic non-parametric transform that makes the 


histogram of T (|Vg\) as similar as hy,. 


Table 2 Algorithm 1 specification 


Algorithm 1: Iterative Specification of the improved adaptive 
regularization method with GHP constraint 


1. Initialization: 
Set the initial estimate as ø” (k) = y , k = 0; 


Set initial regularization parameters a, 27 and step length u. 


2. Loop: iterate on k = 0, 1,..., maxiter 


(a) Let g(k)= Ho’ (k) 

(b) Update G (Vg (k)): 
G(Vg(k)) =sen (Ve (k)):T((Ve (K) 

(c) Define f =G(Vg(k)) 

(d) Compute VJ 


o’(k 
VJ =2H" (y—H"o°(k))+ 2a ZMH) , 29H'V" (f -Vg) 


(e) Update o° 


o'(k+1)=0°(k)+u-VI 
© ke k+1 


until stop criteria satisfied 


There are two regularization parameters (a and 8) in our proposed algorithm. In 
the simulation, a is first determined by Morozov method, % is set empirically. Since 
the parameter 9 balances the regularization term with the histogram preservation term, 
a large 9 is preferred. Actually, % cannot be set too large because of the existence of 
estimation errors in reference histogram estimation. In the simulation, 9 and u are set 
to 1 and 0.25, respectively. For the gradient histogram estimation for Ho, we use 
the same strategy as [23] where the range of f is between 0.001 and 5, and of y is 


from 0.1 to 2.5. 


IV. ALGORITHM PERFORMANCE 

This section investigates the algorithm performance, and focuses on the two crucial 
aspects, i.e., the quantitative evaluation of resolution enhancement and algorithm 
processing accuracy. 

A. RESOLUTION ENHANCEMENT 

Whereas various definitions of effective resolution exist, there are two common 
working definition: the resolving capability for two closely spaced objects and the 
3dB beamwidth of spatial response of point target. The objects are considered 
individually “resolved” if there is a 3 dB change of value between them against a high 
contrast background. To assess the resolution enhancement of the proposed method, 


we take a synthetic o’ field including these features as an example, which is shown in 


Fig. 5(a). In the synthetic o° field, the left squares are 1kmx1km in size and spaced 
2km apart, the right square is 2kmX2km, and the width of the narrow and broad lines 
are Ikm and 2km, respectively. 

The synthetic o? field reconstruction is conducted using Monte Carlo simulation. 
In our simulation, the ‘slice’ measurements of the scatterometer are used to verify the 
resolution enhancement of reconstruction. Based on the measurement model 
described in section II, the measurement result is shown in Fig. 5(b). Fig. 5(c)-(h) 
illustrate the reconstructed o° using three classical scatterometer image reconstruction 
methods (SIR, AART and MART), Tikhonov regularization method (TR) and two Lı 
norm based regularization methods (TV and Lı). The reconstructed result of the 
proposed method is displayed in Fig. 5(i). Noted that due to the equivalence of the 
MAP method proposed by Willams [5] and TR method, the MAP result is not 
descripted repeatedly. 

As expected, the SIR result is excessively smooth and seriously spread out, and 
the noise amplification is least. The AART, TR and TV methods produce better 
resolution enhancement than SIR, but make the output o field become nosier. The Li 
and our proposed method obtain the best resolution enhancement, and our method has 


less noise amplification. 


(a) Truth (b) Measured (c) SIR 


(d) AART 


Fig. 5 (a) Synthetic o’ field, (b) noisy measurement, and reconstructed results of 
(c) SIR, (d) AART, (e) MART, (f) TR, (g) TV, (h) Lı, and (i) the proposed 
reconstruction methods with K,=0.08. 

To future investigate the resolution quantitatively, the profiles of the horizontal 
(azimuth direction) cut of o? values extracted from Fig. 5 across the top row of spots 
are analyzed as shown in Fig. 6. For clarification, Fig. 6(a) shows a comparison 
between the proposed method with three scatterometer image reconstruction methods 
(SIR, AART and MART), Fig. 6(b) presents the comparison between our method 
with the three regularization methods (TR, TV and Lı). From Fig 6(a), it is apparent 
that, except for the SIR method, the other methods realize the resolution of 2km. Our 


method recovers the truth signal with the highest accuracy. We note that, in Fig. 6(b), 
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all methods achieve 2km resolution, and our proposed and L; method have better 
amplitude recovery than the other two. Overall, the proposed method shows good 


performance in resolution enhancement than other methods except for L; method. 
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Fig. 6. Horizontal cuts of o’ values extracted from top row of Fig. 5 


B. DFPSCAT o’ RECONSTRUCTION FRAMEWORK 

To illustrate the application of the proposed deconvolution method for 
DFPSCAT measurements of SWE, simulation is initially used. Fig. 7 is the simulation 
flow chart of DFPSCAT o° reconstruction. In the process, DMRT-QCA model [24] is 
adopted as the snow forward model to obtain the surface original backscattering 
coefficient a, . The DMRT-QCA model takes incidence angle, frequency, 
polarization, snow parameters (density, depth and grain size) and surface soil 
parameters (here it is assumed the surface is frozen ground) as input parameters. The 
observation system outputs the measured backscattering coefficient o based on the 
echo signal model described in section II. Then, the reconstructed backscattering 


coefficient o’ will be obtained after some super-resolution algorithms. 


DMRT-QCA 
Model 


f | 
o°_!| | DFPSCAT | lo? 
System 
| 


| Antenna Noise | 


Fig. 7 Simulation flow chart of DFPSCAT o” reconstruction 


The existing SWE spatial distribution data includes AMSR-E/Aqua L3, ESA 
Globe Snow products, Snow Data Assimilation System (SNODAS) Data Products [25] 
and some airborne and snow pit database. Since AMSR-E/Aqua L3 and ESA SWE 
products have low spatial resolution of 25km, and observations of airborne and 
ground pits have limited spatial distribution, SNODAS SWE dataset becomes the 
optimal alternative for our application. The SNODAS SWE dataset assimilates 


satellite-derived, airborne, and ground-based observations of snow covered area and 


snow water equivalent, and outputs SWE gridded data for the continental United 
States with 1km spatial resolution and 24 hours temporal resolution [25]. 

In the paper, the Montana state of the United States is chosen as the study 
domain to investigate algorithm performance of different SWE distribution. In the 
study domain, the western part is flat (referred to as Coastal Plain), the middle part 
includes gentle hills and valleys (Foothill region), whereas the southern part is 
characterized by mountain ridges (Mountain region). We choose two days SWE data 
(February 28 in 2014 and 2015, respectively) as the input SWE distribution. Since 
actual SWE values are less than 700mm, SWE values above 700 mm are removed. 


The SWE distributions are shown in Fig. 8. 


SWE 20140228 SWE 20150228 mm 


44° ys : | 44° yee 
116° W 114° W112 W110 W108 W106. W104 116° W 114° W 142° w 110 W 108° W 106 W 104 W 


Fig. 8 SWE distributions of the Montana State 
C. ALGORITHM COMPARISON 
To investigate the performance of the proposed method, an algorithm comparison 
among the proposed method and some existing resolution enhancement methods is 
conducted. Based on the former DFPSCAT o reconstruction framework, we take the 
Ku-band VV polarization øo? in Feb. 28, 2015 as a test image, in which three 


different areas, corresponding to the Coastal Plain, Foothill region and Mountain 
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region of Montana State, have been chosen as study regions. Fig. 9 illustrates the test 


image and three study regions marked with red rectangle. 


Test image dB 


Fig.9 Test image and three study regions 
Fig. 10 shows the reconstructed o’ results of the three study regions through 
various reconstruction methods for a particular noise realization of K,=0.08. It can be 
seen that the SIR, AART, MART and TV methods produce smoother image than the 
other three. For TR, Lı and our proposed results, Lı shows the best resolution 
enhancement, but with largest noise amplification; TR has least noise amplification 
whereas sacrifices the resolution; the proposed method shows a comparable resolution 


enhancement as Lı, meanwhile has less noise. 


Truth MART 
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(a) region I 


SIR AART 


L1 


Truth MART 


TR TV 


(b) region II 


SIR AART 
TR TV 


(c) region III 


Truth MART 


Proposed 


Fig. 10 The reconstructed o’ images (in decibels) over three study areas with 
K,=0.08. 
Besides, the quantitative assessment on these super-resolution results is 
conducted. Two evaluation index, correlation coefficient and relative error, are 


computed from the difference between the “truth” and estimated images, where the 


relative error is defined as lo? -0° | /o? 


,/| 


the above methods are shown in Table 3 and Table 4, respectively. 


> .The quantitative experimental results of 


From Table 3, it is obvious that, with the increasing K,, the correlation coefficient 
of reconstructed results gradually decreases. When K,, is less than 5%, the proposed 
method has the highest correlation coefficient. For various K,, the proposed method 
has higher correlation coefficient than L; method, whereas lower than TR which is 
caused by noise amplification. 

Table 4 shows the relative errors results of different methods. It is noted that the 
SIR, AART, MART, TR and TV methods have similar relative errors for the same K,. 
The Lı method has the largest relative errors for different K, values. Compared with 
Lı, the proposed method shows a significant reduction of the error. Moreover, the 


relative error of the proposed method does no increase much with comparison of the 


others. 
Table 3 The correlation coefficient results by different methods 
Region Kp SIR AART MART TR TV Lı Proposed 

0.05 0.85 0.86 0.87 0.88 0.86 0.81 0.89 

0.08 0.84 0.84 0.85 0.86 0.83 0.73 0.84 

I 0.10 0.83 0.83 0.83 0.85 0.82 0.65 0.81 
0.12 0.82 0.82 0.82 0.83 0.81 0.6 0.79 

0.15 0.82 0.81 0.82 0.82 0.8 0.56 0.76 

0.05 0.93 0.93 0.94 0.94 0.93 0.89 0.95 

0.08 0.92 0.92 0.93 0.93 0.92 0.83 0.91 

I 0.10 0.92 0.91 0.92 0.92 0.91 0.78 0.88 
0.12 0.91 0.91 0.91 0.92 0.91 0.73 0.88 

0.15 0.91 0.91 0.91 0.91 0.9 0.7 0.86 

0.05 0.96 0.96 0.96 0.96 0.97 0.82 0.97 

0.08 0.96 0.96 0.96 0.95 0.96 0.76 0.9 

Ill 0.10 0.95 0.95 0.95 0.94 0.95 0.73 0.86 
0.12 0.95 0.95 0.95 0.94 0.95 0.69 0.86 


0.15 0.95 0.94 0.95 0.93 0.94 0.67 0.83 


Table 4 The relative errors (dB) results by different methods 


Region K, 

0.05 

0.08 

I 0.10 
0.12 

0.15 

0.05 

0.08 

II 0.10 
0.12 

0.15 

0.05 

0.08 

Il 0.10 
0.12 

0.15 


SIR 


0.37 
0.38 
0.39 
0.39 
0.4 
0.52 
0.54 
0.55 
0.56 
0.57 
0.2 
0.22 
0.23 
0.25 
0.25 


AART MART 


0.36 
0.38 
0.39 
0.4 
0.4 
0.52 
0.54 
0.56 
0.56 
0.57 
0.21 
0.22 
0.24 
0.25 
0.26 


0.35 
0.37 
0.39 
0.39 
0.4 
0.52 
0.54 
0.55 
0.56 
0.57 
0.2 
0.22 
0.23 
0.24 
0.25 


TR 


0.34 
0.36 
0.37 
0.38 
0.39 
0.48 
0.51 
0.52 
0.54 
0.55 
0.18 
0.22 
0.26 
0.31 
0.36 


TV 


0.36 
0.39 
0.4 
0.4 
0.41 
0.53 
0.55 
0.56 
0.56 
0.57 
0.19 
0.21 
0.23 
0.25 
0.26 


Ly 
0.51 
0.62 
0.71 
0.79 
0.85 
0.67 
0.78 
0.88 
0.96 
1.05 
0.37 
0.42 
0.45 
0.52 
0.56 


Proposed 


0.34 
0.38 
0.43 
0.44 
0.48 
0.45 
0.52 
0.59 
0.59 
0.64 
0.21 
0.27 
0.32 
0.33 
0.36 


D. SWE BACKSCATTERING COEFFICIENT RECONSTRUCTION 


Since the goal of our improved adaptive regularization method is to obtain 


high-resolution o? used for SWE retrieval application, the processing accuracy is 


crucial. In this section, the processing accuracy of the proposed method is investigated 


to validate the feasibility for SWE retrieval. Fig. 11 shows o? images over Montana 


State using one-day data of February 28 in 2015 for DFPSCAT. The figures from left 


to right of each row respectively represent original, observed and reconstructed o? 


images. The figures from top to bottom of each column represent X-band VV-pol, 


X-band HV-pol, Ku-band VV-pol and Ku-band HV-pol, respectively. It is obvious 


that the reconstructed o° images have significant resolution improvement compared 


with the observed ones. To further verify the algorithm performance, another similar 


experiment is conducted using the data of the same day whereas for different years. 
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Fig. 12 exhibits the original, observed and reconstructed o? images of Montana State 


in the year 2014. The same results can be obtained as year 2015. 


X-band VV Original X-band VV Observed X-band VV Reconstruced 


44 N : 44 N3 44 N 
116 W 114 W 112 Ww 110 W 108 W 106 W 104 116 W 114 w 117 w 110 W 108 w 106 W 104\ 116 W 114w 1442 w 110 W 108 w 106 W 104 W 


(a) (b) (c) 


X-band HV Original X-band HV Observed X-band HV Reconstruced 


4 N 
ay ae 114 W 112 W 110 W 108 W 106 W 104 W 
(d) (e) (f) 


Ku-band VV Original Ku-band VV Observed Ku-band VV Reconstruced 


4N 4 N - 4 N= = 
tiw naw 112 w 110 W 108 W 106 W 104 Mw 114 W112 W 110 W 108 W 106 W 104 V “NS w naw 112 w 110 W 108 w 106 W 104 W 
(g) (h) (i) 


Ku-band HV Original Ku-band HV Observed Ku-band HV Reconstruced 


ab Ñ 3 
116 W 114 W 112 w 110 W 108 w 106 W 104 W 


i) (k) 0) 


4 N á 44 N 
“aN W 114 w 112 w 110 W 108 w 106 W 104 116 W 114 W 112 w 110 W 108 w 106 W 104 V 


Fig. 11 The original, observed and reconstructed o° images (in decibels) over 
Montana State using one day data of Feb. 28th in 2015 with K,=0.08. (a)-(c) are 
X-band VV-pol; (d)-(f) are X-band HV-pol; (g)-(i) are Ku-band VV-pol; (k)-(I) 


are Ku-band HV-pol. 
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X-band VV Original X-band VV Observed X-band VV Reconstruced 


44N™ - LNS d d 
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44 e . 
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Ku-band VV Original Ku-band VV Observed Ku-band VV Reconstruced 
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(g) (h) (i) 


4 N 
“aN W 114 W 112 W 110 W 108 W 106 W 104 
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Fig. 12 The original, observed and reconstructed o° images (in decibels) over 
Montana State using one day data of Feb. 28th in 2014 with K,=0.08. (a)-(c) are 
X-band VV-pol; (d)-(f) are X-band HV-pol; (g)-(i) are Ku-band VV-pol; (k)-(I) 
are Ku-band HV-pol 

The relative error and correlation coefficient of measured and reconstructed o° 
images for different K, are shown in Table 5. For each Kp, the first row data is relative 
error in decibels, and the second row data is correlation coefficient .For the day in 
year 2015, the relative errors of the reconstructed o? images are about 0.30dB for 


X-band VV-pol, 0.42dB for X-band HV-pol, 0.27dB for Ku-band HV-pol, 0.33dB for 


Ku-band HV-pol for various K,. They all satisfy the accuracy requirement for SWE 
retrieval except for the X-band HV-pol case. The reason is that the Coastal Plain has 
the highest SWE and relatively larger variation compared with the other regions. The 
corresponding relative errors are about 0.36dB, 0.34dB, 0.22dB, 0.31dB for the year 
2014, and they all satisfy the accuracy requirement. 

In addition, compared with observed results, the relative errors of the 
reconstructed o° images nearly decrease more than 0.03dB for the Ku-band in two 
years. For the X-band, the relative errors decrease more than 0.06dB in two years. The 
correlation coefficient of reconstructed results of both two years shows significant 
improvement than that of observations. 

Based on the above analysis, the proposed reconstructed method could satisfy the 
accuracy requirement of 0.5dB for SWE retrieval and produces high correlated 
reconstructed o° images with the original ones. The results also indicate that the 


proposed method is insensitive to measurement noise and has stable performance. 


Table 5 The relative errors (dB) and correlation coefficient of measured and 
reconstructed o” images 


20150228 
X-band VV X-band HV Ku-band VV Ku-band HV 
O, o; a on O, o; G on 
K.=0.05 0.37 0.23 0.54 0.36 0.28 0.22 0.39 0.27 
ee 0.50 0.86 0.77 0.90 0.87 0.93 0.76 0.88 
K.<0.08 0.37 0.28 0.54 0.40 0.28 0.25 0.39 0.30 
oo 0.49 0.70 0.77 0.86 0.97 0.90 0.75 0.84 
37 : 54 0.42 0.29 0.26 0.39 0.33 
K, =0.10 0.3 0.30 0 


0.37 0.31 0.54 0.45 0.29 0.27 0.39 0.35 


KS 0.47 0.61 0.76 0.82 0.86 0.87 0.75 0.79 
K. -0.15 0.37 0.32 0.54 0.47 0.29 0.31 0.40 0.37 
oe 0.46 058 0.76 0.80 0.85 0.86 0.74 0.76 
20140228 
X-band VV X-band HV Ku-band VV Ku-band HV 
On o’ o’ a O m o On o? 
K. =0.05 0.43 0.28 0.42 0.29 0.24 0.19 0.35 0.26 
pe 0.54 0.88 0.80 0.90 0.86 0.92 0.79 0.90 
K. =0.08 0.43 0.34 0.42 0.32 0.24 0.21 0.35 0.29 
a 0.54 0.76 0.80 0.87 0.86 0.88 0.79 0.87 
K. -0.10 0.44 0.36 0.43 0.34 0.24 0.22 0.36 0.31 
a 0.54 0.71 0.80 0.85 0.86 0.87 0.79 0.85 
K. <0.12 0.44 0.37 0.42 0.35 0.24 0.23 0.36 0.32 
ae 0.53 0.67 0.80 0.84 0.84 0.85 0.79 0.83 
K, =0.15 0.44 0.38 0.43 0.37 0.25 0.24 0.36 0.33 


0.52 0.64 0.79 0.82 0.84 0.83 0.78 0.81 


V. REAL DATA VALIDATION 
The previous analysis is based on simulation processing. In this section, we now 
illustrate the feasibility of the deconvolution framework and the effectiveness of the 
deconvolution algorithm with actual scatterometer data. For the existing pencil beam 
scatterometer data like SeaWinds, the azimuth sampling interval is more than ten 
kilometers which is not satisfied the requirement of high azimuth sampling rate in our 
framework. Therefore, the fan beam scatterometer data become the optimal 


alternative. 


An ASCAT (Advanced Scatterometer on Metop) is a C-band fan beam 
scatterometer, which has two swaths with three beams on each side, as shown in Fig. 
13. Each beam is splited into 256 nodes through range processing, and the 
corresponding measurement values (radar backscatter, incidence angle, azimuth angle) 
along with location (in latitude and longitude) is reported in the “full resolution” 
product (SZF) [26]. We take the right middle beam measurement as an example, and 
combine its multiple measurements for azimuth resolution enhancement. In this case, 
the measurements are 5.6km apart, and the “slice” measurements have an aperture 


function with an effective size of approximately 6kmxX35km. 


ground track 


3 | 
'—550km—— «+ — 360km—»l<——550km—— 
l l l l 


Fig. 13 The ASCAT swath geometry where the beam numbering matches the 


convention in the SZF L1B files. 
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Fig. 14 Images generated from real ASCAT scatterometer data. (a) observation 
image, (b) the reconstructed image processed by our proposed method. 

The costal line of Mogadishu centered at 45°E, 1°N is chosen as the study region. 
Both observation and reconstructed o° images created from the scatterometer 
measurements are shown in Fig. 14. Through the comparison of the observed and 
reconstructed images, it is apparent that the reconstructed image discriminates the 
coastal line clearly and preserves more detail information. It should be noted that 
noise enhancement in the reconstructed image is not obvious, which is mainly due to 
the lower noise level of the ASCAT measurements. 

VI. CONCLUSIONS 
This paper introduces an improved adaptive regularization deconvolution method 
with gradient histogram preservation constraint to enhance azimuth resolution, and 
successfully applies it to the DFPSCAT forward-looking case. Compared with the 
over smoothness of adaptive regularization method, the proposed method can recover 


more details. To illustrate the feasibility of the deconvolution framework and the 


effectiveness of the proposed algorithm, both simulation and real scatterometer data 
processing are implemented. 

In the simulation, a synthetic o’ field reconstruction and SWE o° reconstruction 
are carried out to evaluate resolution enhancement and algorithm processing accuracy 
quantitatively. Meanwhile, the proposed method is compared with some existing 
resolution enhancement methods. The results show that the proposed method can 
obtain a comparable resolution enhancement as L; method. From the quantitative 
analysis of processing accuracy, the proposed method has highest correlation 
coefficient and least relative error when K, is less than 5%, and shows higher 
accuracy than Lı method with higher correlation coefficient and lower relative errors. 
Take SWE distributions over Montana State in two different days as examples to 
judge effectiveness of the proposed method for SWE retrieval application. It turns that 
our proposed method satisfies the accuracy requirement of 0.5dB for SWE retrieval 
and produces the high correlated reconstructed o° images with the truth ones. 
Compared with the observation results, the relative errors of the reconstructed oP 
images nearly decrease more than 0.03dB for the Ku-band in two polarizations, and 
0.06dB for the X-band. 

In addition, we take ASCAT “full resolution” product (SZF) for further validation. 


Compared with observation, the reconstructed image discriminates the coastal line 


clearly and shows better performance. 


ACKNOWLEDGMENT 


The work of this paper is supported by the Strategic Priority Program on Space 


Science from Chinese Academy of Sciences, under the project of intensive study of 


candidate mission of Water Cycle Observation Mission (WCOM) (Grant Number: 


XDA04061202). 


[1] 


[2] 


[3] 


[4] 


[5] 


REFERENCES 
M. S. Lin, Y. G. Zhang, and X. Z. Yuan, "The development course and trend 
of ocean remote sensing satellite," Haiyang Xuebao, vol. 37, pp. 1-10, Jan. 
2015. 
M. W. Spencer, T. Wu-Yang, and D. G. Long, "High-resolution measurements 
with a spaceborne pencil-beam scatterometer using combined range/Doppler 
discrimination techniques," Geoscience and Remote Sensing, IEEE 
Transactions on, vol. 41, pp. 567-581, 2003. 
D. S. Early and D. G. Long, "Image reconstruction and enhanced resolution 
imaging from irregular samples," IEEE Transactions on Geoscience and 
Remote Sensing, vol. 39, pp. 291-302, Feb. 2001. 
D. G. Long, P. J. Hardin, and P. T. Whiting, "Resolution enhancement of 
spaceborne scatterometer data," Geoscience and Remote Sensing, IEEE 
Transactions on, vol. 31, pp. 700-715, 1993. 
B. A. Williams and D. G. Long, "Reconstruction From Aperture-Filtered 
Samples With Application to Scatterometer Image Reconstruction," 
Geoscience and Remote Sensing, IEEE Transactions on, vol. 49, pp. 


1663-1676, May 2011. 


[6] 


[7] 


[8] 


[9] 


[10] 


[11] 


J. S. Shi, X. L. Dong, T. J. Zhao, J. Y. Du, L. M. Jiang, Y. Du, et al., "WCOM: 
The science scenario and objectives of a global water cycle observation 
mission," in [EEE International Geoscience and Remote Sensing Symposium 
(IGARSS) 2014, pp. 3646-3649. 

D. Xiaolong, H. Liu, Z. Z. Wang, J. C. Shi, and T. J. Zhao, "WCOM: The 
mission concept and payloads of a global water cycle observation mission," in 
IEEE International Geoscience and Remote Sensing Symposium (IGARSS), 
2014, pp. 3338-3341. 

P. Clissold, Coreh2o-cold regions hydrology high resolution observation. The 
Netherlands: ESA Communication Product Office, 2008. 

M. W. Spencer, W. Chialin, and D. G. Long, "Improved resolution backscatter 
measurements with the SeaWinds pencil-beam scatterometer," Geoscience 
and Remote Sensing, IEEE Transactions on, vol. 38, pp. 89-104, Jan. 2000. 

J. C. Guan, J. Y. Yang, Y. L. Huang, and W. C. Li, "Maximum a 
posteriori-based angular superresolution for scanning radar imaging," 
Aerospace and Electronic Systems, IEEE Transactions on, vol. 50, pp. 
2389-2398, July 2014. 

R. Sethmann, B. A. Burns, and G. C. Heygster, "Spatial resolution 
improvement of SSM/I data with image restoration techniques," Geoscience 
and Remote Sensing, IEEE Transactions on, vol. 32, pp. 1144-1151, Nov. 


1994. 


[12] 


[13] 


[14] 


[15] 


[16] 


[17] 


[18] 


[19] 


" 


M. A. Richards, "Iterative noncoherent angular superresolution in 
Proceedings of the 1988 IEEE National Radar Conference, 1988, pp. 100-105. 
J. W. Zou, M. B. Zhu, X. P. Li, and W. Dong, "Norm regularization method 
and its application in radar azimuth super-resolution," in TENCON 2013 - 
2013 IEEE Region 10 Conference (31194), 2013, pp. 1-4. 

W. Qian, Y. Ting, and D. Xiaolong, "Super-resolution scatterometer image 
reconstruction using total variation regularization method," in Geoscience and 
Remote Sensing Symposium (IGARSS), 2014 IEEE International, 2014, pp. 
3065-3068. 

M. Y. Zou and U. Rolf, "Circulant and aperiodic models of deconvolution: A 
comparison," Signal Processing, vol. 64, pp. 185-192, Jan 1998. 

A. M. Tekalp and M. I. Sezan, "QUANTITATIVE-ANALYSIS OF 
ARTIFACTS IN LINEAR SPACE-INVARIANT IMAGE-RESTORATION," 
Multidimensional Systems and Signal Processing, vol. 1, pp. 143-177, Jun 
1990. 

M. Y. Zou, Deconvolution and Signal Recovery. Beijing: National Defense 
Industry Press, 2001. 

G. A. Ryzhikov and M. S. Biryulina, "Sharp deconvolution in elimination of 
multiples," Geophysics, 1998. 

W. M. Zuo, L. Zhang, C. W. Song, D. Zhan, and H. J. Gao, "Gradient 


Histogram Estimation and Preservation for Texture Enhanced Image 


[20] 


[21] 


[22] 


[23] 


[24] 


[25] 


Denoising," Image Processing, IEEE Transactions on, vol. 23, pp. 2459-2472, 
June 2014. 

J. K. Patel and C. B. Read, Handbook of the normal distribution. New Y ork: 
New York: Marcel Dekker, 1982. 

W. Guang Xin, L. Hong, and Z. Zhen Rong, "Regularized resolution 
enhancement of point-based features of synthetic aperture radar image using 
variable quasi-norm," in Signal Processing Systems (ICSPS), 2010 2nd 
International Conference on, 2010, pp. V1-500-V 1-504. 

R. C. Gonzalez and R. E. Woods, Digital Image Processing. Beijing: 
Publishing House of Electronics Industry, 2005. 

Z. Wangmeng, Z. Lei, S. Chunwei, and D. Zhang, "Texture Enhanced Image 
Denoising via Gradient Histogram Preservation," in Computer Vision and 
Pattern Recognition (CVPR), 2013 IEEE Conference on, 2013, pp. 1203-1210. 
T. Leung, P. Jin, L. Ding, L. Zhongxin, D. W. Cline, and T. Yunhua, 
"Modeling Active Microwave Remote Sensing of Snow Using Dense Media 
Radiative Transfer (DMRT) Theory With Multiple-Scattering Effects," 
Geoscience and Remote Sensing, IEEE Transactions on, vol. 45, pp. 990-1004 
Apr. 2007. 

N. O. H. R. S. Center, Snow Data Assimilation System (SNODAS) Data 
Products at NSIDC. Boulder, Colorado USA: National Snow and Ice Data 


Center, 2015. 


2 


201703.00309v1 


chinaXiv 


[26] 


R. D. Lindsley, "Estimating the ASCAT Spatial Response Function endnote," 
Brigham Young University, Provo, Technical Report # MERS 14-01, May 


2014. 


